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ABSTRACT 

Rolling element bearings are used in most electrical machines, especially for small 
and medium size applications. Under non-ideal operating conditions, ball bearing 
condition degrades by fatigue, ambient vibration, misalignment, overloading, 
contamination, corrosion from water or chemicals, improper lubrication, shaft currents 
and residual stress left from the bearing manufacturing process. All of these conditions 
eventually lead to increased vibration and acoustic noise during machine operation which 
at some point in time results in unexpected bearing failure. Over the years, a great number 
of publications have been devoted to the detection of mechanical faults, including rolling 
element bearing defects and torsional defects, in electrical machines based on Electrical 
Signature Analysis (ESA). It has been observed that these faults can affect either the stator 
to rotor air-gap distribution or the running speed of the machine, which can be reflected 
in the signature of the electrical signals. However, the physical link between the 
mechanical degradation and the electrical signature is still not explained well. 

A multi-physics model is developed by joining the detailed mechanical model of 
a rotor bearing system and the electrical model of a synchronous machine in this research. 
This combined model is capable of describing the transmission of information originating 
from bearing faults and their impact on the variations of the measured electrical signals. 
The electrical machine model is developed based on winding function approach and its 
validity is demonstrated by a more accurate Finite Element Method (FEM) model. The 


mechanical model consists of a high fidelity rotor-bearing system with detailed nonlinear 


ball bearing model and a flexible finite element shaft model. It is validated using the 
housing vibration data collected from some experiments. 

Generalized roughness bearing anomalies are linked to load torque ripples and 
airgap variations, while being related to current signature by phase and amplitude 
modulation. Considering that the induced characteristic signatures are usually subtle 
broadband changes in the current spectra, these signatures are easily affected by input 
power quality variations, machine manufacturing imperfections and environmental noise. 

In this research, a new algorithm is proposed to isolate the influence of the external 
disturbances of power quality, machine manufacturing imperfections and environmental 
noise, and to improve the effectiveness of applying the ESA for generalized roughness 
bearing defects. The results show that the proposed method is effective in analyzing the 
generalized roughness bearing anomaly in synchronous machines. Furthermore, the 
electrical signatures are analyzed in a synchronous machine with bearing defects. The 
proposed fault detection method employs a Zoomed Fast Fourier Transform (ZFFT) and 
Principal Component Analysis (PCA) and it is also tested on the available experimental 
data. The results show that amplitude induced electrical harmonics are related to the level 
of vibration, and the electrical signatures are affected heavily by other variables, such as 
power quality and load fluctuation. The proposed method is shown to be effective on 


detecting generalized roughness bearing defects in synchronous machines. 
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CHAPTER I 
INTRODUCTION AND LITERATURE REVIEW 
1.1 Overview 

Ranging from high voltage generators used in decentralized power plants including 
electric utilities and district heating plants to low voltage generators which are ideal for 
standalone and emergency power supplies for schools, hospitals, offices and factories, 
synchronous AC generators are widely used in the power generation industry. For 
example, the capacity of installed diesel generators in the United States, reportedly as high 
as 350,000 units with more than 127 gigawatts [1], constitutes a significant component of 
the nation’s power generating infrastructure. 

The reliable operation of power generation systems is very critical because of the 
enormous financial losses resulting from disruption of production, machine replacement, 
human injury and/or fatality, etc. [2]. Although synchronous AC generators are generally 
well manufactured and robust, mechanical anomalies, such as the sub-synchronous 
resonance induced by load troubles, high vibrations resulting from worn bearings, and 
greater torque variation and oscillation caused by faulty prime mover are likely to occur 
in abnormal conditions. These incipient faults, which degrade system performance can 
also lead to failures of other components, finally resulting in an unexpected outage event. 
System operators use condition monitoring (CM), which is a major component of 
predictive maintenance, to monitor machinery and schedule appropriate maintenance 
based on the information that CM provides. By using CM, maintenance planners can 


arrange service during a planned downtime and avoid the costs of an unexpected outage. 


Vibration analysis is the most commonly used technique for CM [3]. Housing 
vibration is measured through a seismic or piezo-electric accelerometer mounted on the 
machine housing. For machines with more critical reliability requirements, radial vibration 
of the shaft is directly measured by an eddy-current proximity sensor mounted on the 
bearing casing. Torsional vibration is measured by a strain gauge directly on the surface 
of the shaft or an encoder attached to it. In a common vibration analysis procedure, the 
vibration data collected and the historical baseline data are compared to observe any trends 
indicative of possible failure. 

Although vibration analysis is quite sensitive to the mechanical anomalies and has 
proved its effectiveness in many situations [4] [5], it can be expensive, difficult to use and 
hard to interpret. While large industrial machines equipped with journal bearings and 
built-in diagnosis sensors can easily be adapted with the technique mentioned above, 
smaller electrical machines that are often equipped with rolling-element bearings and have 
no built-in vibration sensors, require additional sensors. Since the associated accessories 
result in higher production and maintenance expenses for smaller machines, the 
application of CM is limited to large and expensive machines [6]. Furthermore, the 
characteristics of the vibration signal obtained are always related with the sensor’s 
mounting location and orientation. A vibration in a particular part of the mechanical 
component is altered while transmitting through the material, which will separate the 
source from the accelerometer. The requirement of professional knowledge of the 
machine’s vibration transmission characteristics makes it difficult to interpret the vibration 


signal. For example, background noise caused by external devices can interfere with the 
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vibration signal. Therefore, the industry increasingly demands low-cost, non-intrusive, 
and simple-to-use CM technologies and products that can be integrated with the existing 
systems and maintenance planning. 
Electrical signature analysis (ESA) was first developed and patented in 1990 by 

Oak Ridge National Laboratory to overcome the shortcomings of the traditional vibration 
method [7]. ESA is based on the concept that a system’s electrical machine generating 
power/torque also acts as a built-in transducer. Changes in the mechanical/electrical 
system will induce variations in electrical signals, e.g., voltages and currents. These 
signals, although small compared to the total energy provided by the machine, can be 
extracted through modern signal processing technology to characterize the system’s 
overall health condition. ESA is an attractive technology for the following reasons: 

e Equipment-mounted sensors are unnecessary 

e Monitoring a system’s health condition can be done remotely and non-intrusively 

e Continuous monitoring capability for online condition assessment is provided 

e Applicability to both large and small electrical machines, either AC or DC 

e High sensitivity to a range of mechanical disorders, e.g., misalignment, worn 


bearings and gears, loose belts, and degraded power systems 
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Figure 1 Surveys of electric machine failures 


About 60% of the common failures in the components of an electrical machine can 
be attributed to malfunctions of mechanical components, such as worn bearings, cracked 
shafts, etc. Figure 1 reports the results of the IEEE-IAS and ERPI surveys in 1985 and the 
PCIC Survey in 1994 [8]. It is clear that about half of failures are initiated in the bearings. 
Single point/localized defect, and generalized roughness/distributed defect are two major 
classes of faults in rolling element bearings used by most electrical machines. The former 
usually occurs as pits, spalls, or cracks on a raceway or rolling elements, and the latter 
usually occurs because of contamination, lack or loss of lubrication, corrosion, etc., 
Unlike the characteristic defect frequencies of signal point faults, the induced 
characteristic signature is usually a subtle broadband change in the spectrum [9] [10]. 
More studies are needed to understand how this type of bearing fault affects the electrical 


signature. 


Other than bearing faults, torsional vibration is another common issue in power 
generation system. Excessive torsional oscillation of a rotor can lead to fatigue, gear 
damage, power fluctuations, etc. Generally, calculations of the torsional natural 
frequencies occur during the design phase, rather than in the testing phase. As mentioned, 
using a strain gauge or an encoder to monitor torsional oscillation increase the complexity 
of the system and makes it less robust. Hence, cost-effective ESA technology has a great 


potential in this area. 


1.2 Literature Review 

This subsection summarizes the important mechanical failures in electrical 
machines and the corresponding diagnosis techniques proposed in the literature. Modeling 
methods for both electromagnetic circuit and bearing faults are reviewed to give a clearer 
understanding of ESA. 
1.2.1 Diagnosis of Mechanical Faults Using ESA 
1.2.1.1 Bearing Faults 

Under non-ideal operational conditions, rolling element bearings are gradually 
deteriorated by fatigue, ambient vibration with non-uniform forces on the balls, 
misalignment, overloading, contamination, water or chemical corrosion, incorrect 
lubrication, and shaft current and residual stress left from the manufacturing process, 
which eventually lead to increased vibration and acoustic noise during machine operation 
which at some level results in unexpected system failures. According to [11], bearing 


faults can be categorized into single point defect and generalized roughness defect, or into 


localized fault and extended/distributed fault [5, 12]. Limited literature has been reported 
on the subject of diagnosis rolling element bearing faults in synchronous machines using 
ESA, while there are plenty of published papers dealing with bearing fault for induction 
machines in the last decade. Considering the similarities in electrical machines, some of 
the proposed ESA techniques for induction machines can serve as a guideline for this 


research on synchronous generators. 
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Figure 2 Typical vibration signals from localized faults in rolling element bearings [5] 


A single point defect produces sharp impulses when the rolling element passes the 
defected area on raceways. The impulses generate mechanical vibrations at the rotational 
speeds of the defected component. Figure 2 shows the structure and typical vibration 


signals from localized faults in different components of a rolling element bearing, where 
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BPFO, BPFI, FTF, and BSF are the ball passing frequency for outer race, ball passing 
frequency for inner race, fundamental train frequency/cage frequency, and ball spin 
frequency, respectively. Periodically, the vibrations can cause air-gap variations between 
the stator and rotor, which can induce harmonics in the air-gap permeance and flux density 
and result in a fault-related harmonic in the machine’s electrical signatures. In [13], the 


relationship of the bearing vibration to the stator current spectrum is given by 
Sorg =| f tkf, | (1.2) 
where 


f, =(N,/2)f,[1- b, cos(£)/d,] BPFO 
f, =(N,/2)f,[1+ 5b, cos(B)/d,] BPFI 
f, =(N,/2)f, [1 - (b; cos(B)/d,)] FIF 


N, is the number of balls, b, is the ball diameter, d, is the pitch, and 2 is the contact 


angle of the ball with race. To simulate the outer race fault, a hole is drilled on the outer 
race of the shaft end bearing of an induction motor. While this bearing fault is an artificial 
failure which does not occur in realistic cases, the bearing frequencies produced are easily 
found in both the vibration and current spectrums. Many studies have noted the potential 
of utilizing electrical signals to diagnose bearing faults. In [14], an amplitude modulation 
(AM) detector based on high order spectrum has been developed to identify the 
characteristic frequencies of outer race bearing faults. A novel method based on a fuzzy 
logic algorithm for incipient bearing defect detection has been presented in [15]. Wavelet 


packet transform (WPT) and principal component analysis (PCA) are used in pre- 
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processing the signals, and then fuzzy logic is used to distinguish the faults. Electro- 
discharge machining (EDM) has been used to generate localized bearing faults with fault 
diameters of 7, 14, and 21 mils at the inner raceway, rolling element and outer raceway, 
respectively. Conventional frequency analysis methods, however, do not work very well 
with low power machines or under non-stationary loads. To address this issue, Short Time 
Fourier Transform (STFT) and Gabor Spectrogram have been used to detect damaged 
bearings on a permanent magnet synchronous machine (PMSM) [16], although the type 
of bearing faults is not specified. Although the above-stated methods have been verified 
by experimental results, the simulated faults were artificially fabricated offline and re- 
assembled into the motor, which would change the system model and made the results 
doubtful. 

External induced vibration is common in many industrial applications, while 
traditional developed electromechanical fault diagnosis methods in electric machines are 
based on the assumption of a constant radial load on the bearing. Recently [17], has 
addressed bearing fault detection according to the air gap variation model on an induction 
motor under this situation. Both vibration analysis and ESA were applied with the healthy 
bearing and a bearing with spheres removed. The author concluded that the current signal 
is a more robust indicator compared to the vibration signal, because its modulations are 
not affected by external excitation. 

Harmonics of specific frequencies are expected to appear in the current spectrum, 
as a single point defect occurs in the bearing. However, when the small localized spalls 


have become extended and smoothed by wear, no sharp impacts will be generated 


necessarily [5]. Instead, a broadband change in the current spectrum will be found. Several 
papers [18, 19, 20, 21] investigated applying ESA on the detection of this generalized 
roughness bearing defect. The authors simulated in-situ bearing deterioration by injecting 
high AC current externally. A model-based technique for bearing fault detection of 
induction motor has been proposed in [19]. The harmonics unrelated to bearing faults were 
filtered out and the residual was estimated via an autoregressive (AR) model. The 
parameters of the AR model increased as the bearing degraded. However, the proposed 
indicator showed poor performance with large load variations. In [18], a Wiener filter 
based noise cancellation technique has been employed to reduce the significant harmonics 
that are unrelated to bearing faults. The author in [20] extended this method by adapting 
it to variable load and running speed. This method has been further proved proficient in a 
claw-pole synchronous generator [21]. The results of vibration or current signatures under 
different external radial loads applied to the shaft investigated in [22], have shown that 
radial load has a negligible effect on the output torque, vibration, and current signals. 

Several studies have focused on the detection and diagnosis of rolling element 
bearing faults in electrical machines based on ESA. However, the physical links between 
the mechanical vibrations and electrical spectrum components is still unclear [23]. Thus, 
there is a strong need to develop an electromechanical coupled model that can explain the 
physical link between vibrations and electrical signatures. 
1.2.1.2 Torsional Vibration 

Torque dynamic variations, torsional oscillations, and load troubles induced by 


internal or external factors are common factors that can be encountered during the 


operation of electromechanical system. Under these conditions, premature failures will 
occur when the cumulative stress acting upon the mechanical components exceeds design 
limitations. For instance, in an engine-generator system, the reciprocating motion of the 
piston attached to the crankshaft mainly generates the torsional vibration. Should a misfire 
occur, which causes an imbalance in the cylinder pressures, the reciprocating frequency 
of the malfunctioning cylinder induces torsional vibration. Although a torque meter or 
strain gauge [24] is generally used to investigate the mechanical behavior of the torsional 
vibration, electrical signature analysis is also considered as a good alternative for 
mechanical monitoring. 

Detecting shaft torsional vibration using stator currents has been first published in 
[25]. This paper has claimed that the current components at specific frequencies which are 
the modulation products of the power system frequency and the vibration frequency were 
the proper indicators of torsional vibration. A general closed form proof and experimental 
validation have been given for both induction and synchronous machines. In [26], 
harmonic analysis of stator current for Permanent Magnet Synchronous Machines 
(PMSM) has been presented to develop a robust real-time rotor unbalance estimation 
scheme for condition monitoring. Forced steady state response of torsional vibration of 
induction machine has been analyzed under different excitation frequencies and 
magnitudes [27]. In [28], a 2DOF dynamic model of a single stage gearbox, which 
accounts for the characterization of the pinion-wheel contact point with damping and 
stiffness, has been used in the electromechanical model to show the effects of external 


input and mesh frequency on the electromagnetic torque and stator current signature in the 
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induction machine. The results showed that both the current spectrum and the estimated 
electromagnetic torque gave significant information about system health . 
1.2.2 Electromagnetic Circuit Modeling 

In an electrical machine operating under ideal conditions, the air-gap between the 
stator and rotor is uniform and the radial magnetic forces are well balanced. Nonetheless, 
due to shaft deflection, inaccurate mounting of the rotor with respect to the stator, bearing 
wear, stator core looseness, etc., the spinning center of the shaft may not be concentric 
with rotor axis and stator bore axis [29]. This will disrupt the fine balance between the 
magnetic forces of adjacent poles, by increasing the unbalanced magnetic pull (UMP) and 
placing more load on the bearing. Excessive stresses resulting from the eccentricity fault 
can cause rubbing of the stator and rotor, eventually damaging the stator core and stator 
windings. When eccentricity occurs, varying inductances cause unbalanced magnetic flux 
within the air-gap that creates fault harmonics in the line current, which can be identified 
in the spectrum. Because bearing faults will cause radial motion of the shaft which varies 
the air-gap length and generate eccentricity, it is feasible to apply an electromagnetic 


model for eccentricity faults on bearing failures. 
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Figure 3 Static eccentricity 


In general, there are two types of eccentricity fault: static eccentricity and dynamic 
eccentricity. Figure 3 shows the case of static eccentricity, where the rotation center is 
coincident with the rotor center but outside the stator bore axis. The minimum air-gap 
position does not change as the rotor rotates. Figure 4 shows the case of dynamic 
eccentricity, where the rotor spins around its own centerline and also whirl around the 
stator axis. The air-gap length varies as the rotor rotates. In reality, both static and dynamic 


eccentricities occur simultaneously, which is known as a mixed eccentricity fault [30]. 
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Figure 4 Dynamic eccentricity 


Many studies have investigated the use of a current signal as the medium to detect 
air-gap eccentricity faults. Some studies have focused on the detection of fault-related 
signatures in the current signals, while others have focused on modeling and simulating 
faults. 

The winding function approach (WF A) is the most widely used method because it 
requires less computation complexity and it is comparably accurate. In [31] [32], the 
sidebands of either the slot frequencies or fundamental have been used to detect air-gap 
eccentricity. However, the underlying theory behind had not yet been established. Toliyat 
[33] was the first to use the winding function approach (WFA) to study the performance 
of induction motor with SE fault. Unlike a traditional lumped model using only the 
fundamentals, WFA accounts for the space harmonics induced by non-sinusoidal 
distributed winding, therefore making it possible to simulate motor performance with SE. 
The simulation results closely agreed with the experiments. 

An improved analytical approach known as the modified winding function 


approach (MWFA) has been proposed in [34] [35]. The simulation as well as experiments 
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of a single phase and a three-phase induction machine under mixed eccentricity faults 
demonstrated the possibility that monitoring slot harmonics in the stator current to detect 
eccentricity faults. MWFA has been further improved by accounting for the skewing of 
the rotor bars and the resulting linear rise in the MMF across the slots [36]. In [37], a 
more practical case, the axial non-uniformity of air-gap known as inclined eccentricity has 
been reported. In [38], startup transient simulations have been performed on an induction 
motor with mixed rotor eccentricity fault using Time-Stepping Finite Element (TSFE) 
approach. In the simulation, the applied voltage was used as an input. The results showed 
that the sidebands of fundamental were well indexed for detecting mixed eccentricity. 

Since the 1990s, the needs to improve energy efficiency and the advancements in 
power electronics have prompted researchers to turn their attention to synchronous 
machines. For example, in [39], a MWFA based model has been created to study DE in a 
salient pole synchronous machine, where the 17" and 19" harmonics have been reported 
as a good indicator of DE fault. By incorporating the open slot effect, a more accurate 
model has been reported in [40]. Rotor faults in turbo-generators have been studied in [41] 
using transient FE analysis. The response of a two pole turbo-generator to DE fault has 
been analyzed in both time and frequency domain. Frequency domain signatures in the 
no-load voltage spectrum of the generator were found to identify the existence of rotor 
fault, even though fault discrimination could not be achieved. 

The different effects of SE and DE on the external electric variables of an on-board 
ship salient-pole synchronous generator have been presented in [42]. Unlike monitoring 


particular fault-related harmonics, the increment of a no-load EMF space vector has been 
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investigated by 3D FEM as a useful tool to detect eccentricity faults. A different oval 
shaped EMF space vector was found with SE but not DE, proving that the shape of the 
loci can be used to distinguish the type of eccentricity. However, no discussion regarding 
the effect of ME was reported. In a later publication [43], compared comparison of an 
MWFA-based model with an improved air-gap function and the 3D FEM model has found 
similar results; the author concluded that no-load voltage harmonics also showed good 
dependence on level and type of eccentricity. In [44], another comparison between MWFA 
and FEM on a 2 kW, 3 phases, 60Hz, star connected salient pole synchronous motor, 
where the fault related harmonics have been found in both models, has been verified by 
experimental results. In [45] and [46], based on a MWFA model, general criteria to predict 
the fault-related voltage harmonics for no-load synchronous generator affected by static 
eccentricity have been carried out, and followed by the use of a 15kVA four pole 
fractional-slot synchronous generator for experimental validation. The SE was attached 
by an innovative flange with the capability of adjusting amplitude and direction of both 
SE and DE. The influences of different excitation schema, AC and DC exciting currents, 
have been discussed. SE has been found to introduce a double fundamental frequency 
ripple in the field current. While the experimental results confirmed the theoretical 
prediction, they showed that linear models useful for the simplest cases were not sufficient 
in providing precise quantitative results when a mixed type eccentricity fault occurs in a 


real machine. 
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1.2.3 Bearing Modeling 

In studying bearing faults and fault mechanisms, the underlying physics and 
behavior of bearing systems should be well interpreted. This understanding is achieved by 
reviewing and developing bearing models, and understanding how experiments can 
validate the modeling methods. 

References [47] [48] have provided a solid foundation for the modeling of bearings 
by developing bearing element stiffness, deflections, damping, constraints, forces, 
moments and general governing equations of motion. All of the derivations were based on 
fundamental theories ranging from Hertzian contact theory to Lagrangian dynamical 
derivations. 

Gupta [49] has provided a detailed formulation of specific bearing dynamical 
behavior which was previously neglected, such as ball-race interactions, ball-cage 
interactions, churning and drag losses, and lubrication effects. In [50], Gupta (1985) 
refined his model and validated it experimentally. Both the simulation and the experiments 
indicated a certain critical shaft speed at which the cage begins to whirl. A more detailed 
analysis has been presented by Datta and Farhang [51], which treated each element in the 
bearing as a nonlinear spring. 

To model the bearing defects, McFadden [52] has taken a different approach by 
using impulse train to model the vibration produced by the impact of a bearing element 
with a single point defect. A multiplying factor on the impulse function was applied to 
represent the level of severity of the defect. Later, the author extended the model to adapt 


for multiple point defects [53]. An experiment was also performed to verify the model. 
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In [12], to study the extended fault located at the inner/outer race of rolling element 
bearings, a combined dynamic model for gears and bearings has been proposed to obtain 
a better understanding of the physical link between the two components. The extended 
faults do not necessarily produce high frequency impact responses, but do modulate the 
gear-mesh signals. Gaussian noise filtered by low pass filter was generated to simulate the 
rough surface. The simulation results mostly agreed well in with the experimental data, 
showing that the extended fault could be clearly detected after removing the gear 
contribution and comparing the cyclic spectra. 

In 2015, Linkai Niu et al. [54] have carried out a systematic study of ball passing 
frequencies (BPFs). Localized surface defects were modeled by considering the ball size, 
the clearance due to material absence, and contact force directions. The study showed that 
several elements, including the shaft speed, external loads, friction coefficient, raceway 
groove curvature, initial contact angle, and defect sizes, affected the BPFs. The 1-2% 
variation of the spacing between impulses caused by a defective ball bearing reported in 
previous studies, was only satisfied when the skidding effect in the bearing is small. Two 
experiments were performed to verify the findings. 

In 2015, an extensive review for the vibration modelling of rolling element 
bearings with localized and extended defects has been undertaken by Sarabjeet et al. [55]. 
The authors also discussed the benefits and limitations of algorithms for estimating the 
size of bearing defects. They recommended conducting a full parametric model including 


load, rotational speed, clearance, and various defect types in the future. 
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1.3 Research Objective 

Over the years, a great number of works have been devoted to detection of 
mechanical faults including rolling element bearing defects and torsional defects in 
electrical machines based on ESA. It has been discovered that these faults can affect either 
stator to rotor air-gap distribution or running speed of the machine, which can be reflected 
in the signature of electrical signal. In addition, many talented works have been 
contributed to the understanding of the dynamics of rolling element bearing while a 
localized or a distributed fault occurs. However, the physic link between the mechanical 
signature and the electrical signature are still not clear enough. The objectives of this 
research can be listed as follows: 

The published literature concerning the use of ESA to detect mechanical faults 
including rolling element bearing defects and torsional defects in electrical machines has 
determined that the faults, which can affect either stator to rotor air-gap distribution or 
running speeds of the machine, can be reflected in the electrical signature. The physical 
link between the mechanical signature and the electrical signature, however, has not been 
determined conclusively. Therefore, the objectives of this research are: 

e To include both rolling element bearing dynamics modeling and synchronous 
generator modeling into a new comprehensive model capable of describing the 
transmission between bearing faults and electrical signal variations 

e To analyze the impact of a synchronous generator with a bearing defect on 
different components in order to guide the development of an ESA-based bearing 


fault detection schema 
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e To conduct an in situ experiment and analyze the electrical signatures caused by 
the general roughness bearing defects to gain practical insights into the use of 


electrical measurements as a bearing fault detection source 


1.3.1 Proposed Approach 

The proposed scope of this dissertation will include synchronous generator 
modeling, bearing modeling, multi-body dynamics modeling for drive trains, and model 
verification. The first step is to model an electrical machine with the capability of 
introducing mechanical defects. As mentioned in the literature review, WFAs are widely 
used in modeling internal faults including short windings, broken rotor bars, eccentricity, 
etc. Compared with the traditional lumped model, WFAs account for all of the space 
harmonics and interpret the influences from the variation of air-gap between stator and 
rotor. A promising synchronous generator model based on MWFA will be built to simulate 
various rotor eccentricity scenarios. A more accurate FEM model created with the 
commercial FE package JMAG will be used to validate the model. 

The majority of the literature on eccentricity faults assumes that the rotor rotates 
along a round circle orbit synchronized with the shaft spinning speed. However, it is not 
the case for realistic machines. As shown in [56], the measured shaft orbits are irregular 
and altered depending on load conditions. To illustrate the practical air-gap variations 
induced by bearing faults, a new rolling element bearing model will include a flexible, 
physically based introduction of mechanical faults into the system. The proposed model 
will provide a quantified system response to a comprehensive set of typical bearing 


failures. 
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The electromagnetic circuit model and the mechanical models including bearing 
defects will be coupled to explore the electrical signatures induced by the faults 
introduced. ESA will be carried out using the simulated data obtained from the inter- 
coupled synchronous machine system and the experimental data from the in-situ bearing 
deterioration experiment. The results will be compared to the trend of vibration signatures. 
1.3.2 Anticipated Contribution 

Reliable operations of synchronous machines are critical in power generation 
systems. Moreover, there is a growing need for inexpensive, non-intrusive, easy-to-use 
diagnostic technique. Although various ESA-based algorithms have been developed and 
verified in the laboratory, a clear theoretical understanding of the physical link between 
the mechanical signature and the electrical signature has not been determined. The 
proposed inter-coupled simulation model will bridge the gap between electromagnetic 
circuit modeling of electrical systems and the dynamic modeling of mechanical system. 
The research will make the following contributions: 

e A new model combining a detailed bearing model in a synchronous generator will 
allow introducing bearing faults 

e Comprehensive analysis of synchronous generator voltages and currents under 
bearing defects that induce air-gap variation 

e A useful tool for analyzing the performance of ESA-based diagnosis algorithms 


in a controlled environment 
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CHAPTER II 
MODELING SALIENT POLE SYNCHRONOUS MACHINES WITH 
ECCENTRICITY FAULTS 

In the last few years, more and more researchers have shown great interests in the 
condition monitoring of synchronous AC machines, on-line and/or off-line, because of 
enormous cost on parts replacement and production lost. If properly used, condition 
monitoring (CM), can extend the lifecycles of various types of ac machines. CM is based 
on modeling the machine faults that are subsequently used to analyze the full effects on a 
machine’s behavior. To analyze the behavior of synchronous machines under faulty 
conditions, developing simple models for the machines which can take the faults into 
account is very important. 

In this chapter, a model based on WFA that can simulate synchronous machines 
under rotor eccentricity faults is built firstly. Secondly, the effects of dynamic and static 
eccentricity are studied using the developed model on the self-magnetizing and mutual 
inductances of a salient pole synchronous machine (SPSM). Thirdly, a finite element 
model built in a commercial software, which is widely employed in the industry, is used 
to validate the WFA model. Lastly, specific harmonic components in the current 


signatures are discussed under different eccentricity level in both WFA and FEM models. 
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2.1 Winding Function Approach 

Conventional modeling methods of electrical machines do not take the winding 
configuration and geometry of the machine into consideration, which results in the 
incapability of simulating the space harmonics changes induced by various fault 
conditions such as broken bars, short turns, and rotor eccentricities. Although the novel 
method FEM provides precise modeling, its application is complicated and the 
computation cost is huge especially for the transient analysis of electrical machines with 
asymmetry in the motor body [57] [58]. WFA is based on the geometry and the winding 
distribution of the machine, which makes it possible to analyze the machine's behavior 
under winding and rotor position variations due to faulty conditions [57]. Because of its 
moderate accuracy and much less computation effort requirement, this method has been 
broadly applied to analyze broken rotor bars [58], fault condition in stator windings [61], 
eccentricities in induction and synchronous machines [35] [44]. In a recent publication 
[62], it was extended for radial and axial non-uniform air gap in SPSMs. 

Calculating the machine inductances accurately is the most important step to get a 
validity model, which can improve the accuracy of the analysis of electrical machines. In 
[63] an analytical equation of the air gap for an eccentric SPSM was presented. The turns 
function and the inverse air-gap function are expanded into Fourier series. The inductances 
are calculated with the simplified analytical solutions. However, the rotor motion is 
assumed to be a synchronous whirling under eccentricity faults. In this research, air gap 


function is not only a function of rotating angle but also a function of rotor offsets. 
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2.1.1 Winding Function 


Figure 5 Winding Function of an elementary SPSM 


The basic idea of winding function approach has been introduced in [64], and the 
derivations below is following the book. An elementary two-pole salient pole synchronous 
machine (SPSM) with four concentrated stator coils A,B,C and D is shown in Figure 5. 
Each stator coil is configured to full pitch. The current which goes into the paper is denoted 
with a cross “x” and is defined as positive. The current which goes out of the paper is 
denoted with a dot “e” and defined as negative. A reference angular position is selected 
as 0 = 0. The path 12 is taken across the gap from stator to rotor at 0 and the path 34 
returns across the gap at 0 + @. 

The same current i is carried in all coils of the windings enclosed by the closed 
path 1-2-3-4-1. According to Ampere’s Law, the following relation can be obtained 


considering the closed path 1-2-3-4-1 illustrated in Figure 5: 
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$H -ai =n(9,0)i, (2. 1) 


12341 
where H is the magnetic field intensity and dl is defined to be along the closed path 1-2- 
3-4-1. This equation states that the circulation of magnetic field intensity in free space 
around any closed path is equal to the total current enclosed by the path. The function 
n(Q, 0) , which expresses the number of turns surrounded by the closed path 1-2-3-4-1, is 
called the turns function. It is a function of ø for stationary coils in the stator. For the 
coils of the rotating rotor, it is a function of g and rotor position 0 . Turns are considered 
as positive, if the current they carry goes into the page. Otherwise turns are considered as 
negative. The turns function n(g,@) for the arbitrary windings is shown in Figure 6, 
assuming that each side of the coil has N turns. 
When the closed path is split into four components, the magnetomotive force 


(MMF) drops in the closed path can be written as: 


MMF,, + MMF,, + MMF, + MMF, = n(9,0)i (2. 2) 
n(¢, 6) 
4N 
3N 
2N 
N 
A B G D A' B' G D 9 


Figure 6 Turns function for the winding distributions 
24 


Because the magnetic permeability of the stator iron is much larger than that of 
the air gap, the iron permeability is assumed to be infinite by ignoring the saturation 


effects. Therefore, the MMF drops MMF,, and MMF, can be neglected and (2. 


2)reduces to: 
MMF,, (0,0) + MMF, (6,9) = n(Q, @)i. (2. 3) 


By Gauss’s law for magnetic fields: 
B-dS=0 
> 2.4 
f (2.4) 


where B represents the magnetic flux density and S represents the surface of an closed 
cylindrical shape. This equation states that there are no magnetic flow sources, and the 
magnetic flux lines always close upon themselves, i.e., the magnetic flux through any 
closed surface is zero [57]. Since all flux is assumed to pass through the air gap and the 


iron permeability is infinite, (2. 4) can be written as: 


2m L 


fB -as = | [utto 0p = syr, | H(@,0\d0 =0 (2. 5) 
S 0 0 0 


where 4 is the free space permeability, rs is the stator inner radius, L is the effective 


axial stack length of the machine. Since B does not vary with respect to the axial length, 


after eliminating the constant components, (2. 6) can be written as: 


2a 
MMF (9,0) _ 


J H(9,0)do = J 20703 Q. 6) 
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Because MMF can be expressed as the product of the magnetizing force and the 


flux radial length, substituting equation (2. 3) into (2. 6) yields: 


ee 0) + MMF (9,0) 1 € ^n(o,0) 


ido 207 
g(9, 0, x, y) o (9,0, x, y) vp 


Substitute the relation of (2. 6) into (2. 7). The equation above can be reduced to: 


í MMF,,(0,0) MMF (9,0) jg- 
o 8(P,9,x,y) — g(g,0, x, y) 


(2. 8) 
MMF, (0.0) Í E yy! - [ 20.8 
o 8(9. 0. x, » 
Therefore, MMF,,(0,@) can be expressed as: 
2x 
n(qo,0 id 
MMF, 0,0) = | — 99 _, ie] | [27(9.0..) do (2.9) 
o 8(9.0. x x») 
By substituting (2. 9) into (2. 3), MMF,,(o,0) can be solved as: 
1 2x 
MMF. (9, 0) = (n(o,0 n(o,0)g (9,0, x, y)do)i 
(9.6) = (n(9. 0) — = falo, 6g (9,0. x. »)do) TN 


[2 (0.6.49 ° 
0 


Therefore, the modified winding function M (9,0) [39] with non-symmetric air- 


gap can be defined by: 
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1 2m 


M(9,0)=n(9,0)- x [n(o.0)g  (9.0.x. y)do 
fe (0.06.x, ydo ° 
0 


(2. 11) 


The general form of the modified winding function is defined in (2. 11). It is 
applicable to any given winding configurations or non-symmetric models. 
2.1.2 Eccentricity Modeling 

Flaws in the manufacturing process, shaft misalignment, unbalanced rotor, and 
vibration excited by connected mechanical components, all these factors may result in 
eccentricity fault in an electrical machine. As shown in Figure 7, there are three types of 
air-gap eccentricity [58], which are static eccentricity (SE), dynamic eccentricity (DE) and 
mixed eccentricity (ME). In the case of SE, the rotor centroid and the center of rotation 
are identical, but they do not coincide with stator centroid. The minimal air-gap location 
is fixed with respect to the stator. Dynamic eccentricity occurs when the rotor centroid is 
off to the center of rotation and the stator centroid coincides with the center of rotation. 
The minimal air-gap location moves while the rotor rotates. When both SE and DE exist, 
ME is generated. The rotor centroid, the center of rotation, and the stator centroid are not 


coincided with each other. 
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Figure 7 The three eccentricities in a synchronous machine. SE, DE, ME 


Apparently, compared to the symmetrical conditions, the winding functions of the 
stator windings and rotor winding of the machine do not make any difference in non- 
uniformity conditions. Nonetheless, the functions of the air-gap length and the mean radius 


will become different in the symmetrical case. In the ME case, the air-gap and mean radius 


functions are illustrated by: 


g«(6,0,x, y) = g($,0)(1— 
(2. 12) 
ôs cos(ġ ~ Ps) 7 ógcos($ i Ps ~ 8)) 
where ge (6, 0) is the airgap function with eccentricity, g (45, 0) is the airgap function for 
a healthy machine, 0 is the rotor position angle in stator reference, $ is the angle in stator 
reference frame, (x, y) are the rotor offsets on vertical and horizontal direction, @, is the 


angle by which the rotation and stator centers are off, and 6, and 6g are the static and 


dynamic relative eccentricity levels which are also functions of (x, y), respectively. 
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2.1.3 Calculation of Inductances 
Calculation of inductances is based upon three assumptions: 
e Slotted stator effects are negligible, 
e Magnetic hysteresis is negligible, and 
e The iron core is not saturated. 


The MMF drop in the iron core is neglectable, because the relative permeability of 
the iron is assumed to be infinity in the model. The MMF distribution can be simply 
calculated using the winding function defined in (2. 11) with the current flowing in the 
coils. 

Let's assume the modified winding function of winding A is M4(, 0, x, y), and 
winding B has a turns function ng(@, 0). According to [64], the flux crossing the air-gap 


at any position by winding A is given by: 


do 


d = MMF (Q, Our, L —— À — 
Z(Q, 0, x, y) 


(2. 13) 


The total flux linkage for winding B due to current in winding A can be defined 


as: 


A 3 frig (@,€F,(0,6)—2? 
= uL ng, (9, 0) F (o, 9) ———— —— 
i ° k=l 0 N " g(9,0, x, y) (2 14) 


2x 
dp 
= AL | n (p, OF (9,0) — — 
i | í C 80x y) 
di 
where n5(9,0) = yon m (9,0) is the turns function for winding B under the assumption 
k=l 


that the g, coils are connected in series. Then the mutual inductance L,, can be calculated 
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by: 


A rL d 
Ly, == |n, OF ,(p,0)—"> — (2. 15) 
Ly Ly 0 g(9, 0, x, y) 


By substituting F(9,0) =M (9,0), into (2. 15), 


A x d 
Dus = -= Url [ns(o.0)M ,(9,0) — — (2. 16) 
4 g(9,0, x, y) 


A 
Using the same procedure, the self-magnetizing inductance of winding A can be 
obtained as: 
2x do 
La = uL | n,(9,0)M (9.0) —À — Q. 17) 
0 g(9,0,x, y) 
Equations (2. 16) and (2. 17) can be utilized to calculate any electrical machine 
mutual and magnetizing inductances even when the machine is experiencing failures such 


0) 


as rotor eccentricity and shorted rotor field. It is because the turns function " (p, and 


M(o,0) 


the modified winding function are calculated in their general forms to include 


even the failures and asymmetries, which can be associated to either stator or rotor. 


2.2 Electrical Circuit Model for Synchronous Machines 
2.2.1 Stator and Field Equations 
All of the inductances for solving the machine equations can be calculated using 
the modified winding function approach (MWFA) described above. Once the inductances 
of the synchronous machine are obtained, the dynamic response of the system can be 
computed. 
The voltage equations for the stator circuits of a typical 3 phase synchronous 


machines can be written as: 
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dA, 


Vanes = Rslapcs + dt (2. 18) 


By the definition of inductance, the stator flux linkages are given by: 
As = Lsiabcs + Lsgiy 
Vapcs and igncs , which are the stator voltage vector and current vector, 


respectively, can be written as: 


Vas las 
Vapcs — | , labcs = le] (2. 19) 


Ves les 


Stator inductance matrix L, is a symmetric 3 by 3 matrix that can be written as: 


Len, y, 8) T ls Lay Qo y, 8) Lac (x, y, 0) 
Ls = Lpa (x, y, 8) Lpp (x, y, 0) + ls Lyc(x, y, 0) (2. 20) 
Leay, 8) Loy (x, y, 0) La y; 0) + ls 


where Li; is the self-inductance of phase i, L;; is the mutual inductance between phase i 
and j, and l, is the leakage inductance of the stator winding. 
Lsf is a 3 by 1 matrix comprised of the mutual inductances between the stator 


windings and the field windings that can be written as: 


Las (x,y, 8) 
Lsf = Lpr (x, y, 8) (2.21) 
Leg (x, y, 0) 


The self-inductances and mutual inductances are calculated by the rotor offset 
(x, y) and the rotating angle 0. 
The stator resistance matrix R, consists of the stator resistance. The resistance in 


each stator winding identical to r;is given by: 
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rn o 0 
Rea r ] (2. 22) 
0 0 5 
The voltage equation of the field winding is given by: 
dy 
(2. 23) 


Ap = Lrslapcs + Leif 
where Rẹ is the field winding resistance, i; is the field current, Lfs = Liy is the matrix of 
the mutual inductances between the field and the stator windings, and Lp is the self- 
inductance of the field winding. 
2.2.2 Mechanical Coupling with the Torque 


Assume the rotor of the synchronous machine is a simple disk with moment of 


inertia J. The equation of motion is given by: 


T, — Tj, =JÖ + BÓ (2. 24) 
where Tę is the electromagnetic torque, T; is the load torque, B is the rotor bearing 
damping, and @ is the rotor position. 

According to the co-energy method, the electromagnetic torque can be calculated 
by 


T 
^ do 


(2.25) 


In a linear magnetic system, the co-energy is equal to the magnetic energy stored 
in all inductances. Thus, when the damping winding is ignored, the co-energy can be 


expressed as: 
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1 
Wo = gU Lsls + If Lelg + If Lgsls + 1g Lsplp) (2. 26) 


Substitute it back to the torque equation, we have: 


aL, 
08 


T, = s(t I, + IF i +217 ^ I) (2.27) 
2.2.3 Assembling the Equations for the Dynamic Model 

To solve for the unknown variables, a set of differential equations, including the 
stator and field voltage equations and the equations of motion for the rotor given above is 
assembled. For the motor, the unknown variables are the currents. When the electrical 
machine runs as a generator, the stator voltages are also unknown. In addition to the 
unknown voltages, the electromagnetic torque becomes the load torque, and the torque 
applied on the shaft becomes the input torque for the mechanical system. The assembled 


set of differential equations is: 


iai 

A7 

do 1 

Gat ze — T;) (2. 28) 
do — 

dt ^ 


A — LI, then the first differential equation becomes: 


d^  dLI UU C — 
dt^ dt dt "dto GM) 
Inductance matrix L is a function of 0, x, y; therefore: 
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dL 00 0L @x@L ðyðL aL 


u au eae decay OB 


Substitute it back to (2. 28), the final expression of the differential equations used 


to find the unknown currents becomes: 


d any gp ge yo yt 
dp t sopra ay) 
GD cu 

dt ; V L) 

do _ 

dp 
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(2. 30) 


(2. 31) 


2.3 Simulation Results and Model Validation 


Figure 8 Cross section model of the SPSM 


To validate the proposed method, a 13kVA, 3 phases, 4 poles, 60Hz SPSM is 
modeled using both WFA and FEM. The main parameters are listed in Table 1. The 
commercial software JMAG is used to create the FEM model. The cross section of the 
modeled machine is shown in Figure 8. The next sections describe a comparison of the 
inductances calculated by both methods followed by a comparison of the voltage/current 


spectrums in both open circuit cases and loaded cases with different dynamic eccentricity. 
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Table 1 Main Parameters of the SPSM 


Name Symbol | Value 
Synchronous speed Ws 1800rpm 
Phase m 3 

Pole pairs D 2 

Rated power Prate | 13kVA 
Number of stator slots N; 36 
Number of slot layers q 2 
Number of parallel circuits C 1 
Number of turns per phase connected in series Nph 124 

Coil span Cspan |6 
Number of turns per pole of field winding N, 76 
Minimum air gap 9o 1.4mm 
Maximum air gap Üm 11.2mm 
Stator stacking length l; 80mm 
Pole arc Parc | 120° 
Stator resistance per phase R; 1.020hm 
Field resistance Ry 0.480hm 


2.3.1 Calculated Inductances 

The accuracy of the electrical circuit model is determined by how well the 
inductances are calculated. Figure 9 depicts the self-inductances and mutual inductances 
over 180 degrees of rotation angle as computed by both WFM and FEM. The rotor is 
moved away from its center 0.7mm (50% of the minimal air-gap) toward the x-axis 
direction. The shapes of Laa, Lap, and Lac are identical, but the phases are shifted by 


1/6. In addition, Laa is levitated by about 0.2H because of the eccentricity along the x- 
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axis. As expected, the shapes of inductances calculated by WFM are smoother than those 


calculated by FEM, because simplifications are made to the magnetic circuit in WFM. 


WFM vs FEM 


Inductance (H) 


200 


Rotor position (degrees) 


Laa-WFM Lab-WFM Lac-WFM Laf-WFM ———Lf-WFM 


Laa-FEM 


Lab-FEM 


Lac-FEM 


Laf-FEM ——=Lf-FEM 


Figure 9 Inductances over different rotating angles 


When the dynamic eccentricity is changed from 0% to 50% for the modeled 
SPSM, the computed inductances are compared as shown in Figure 10. Note that the 
absolute values of the inductances become larger with the increase of the eccentricity, and 
the inductances increase nonlinearly with the eccentricity changes. The inductances’ 
errors between the two modeling methods are shown in Table 2. Because the purpose is 
to investigate the influence of airgap variations on electrical signals, we are trying to find 
a balance point between accuracy and computational cost. According to the comparison, 
the errors of the field-related inductances are much higher than the stator-related 
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inductances, which are still acceptable when they match the curve pattern of the 


inductances in FEM. 


Laa Lab 
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0 0.2 0.4 0.6 
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Figure 10 Inductances vs eccentricity 
Table 2 Errors between the inductances computed by WFM and FEM 
Offset | 0 0.05 | 0.1 0.15 | 0.2 0.25 | 0.3 0.35 | 0.4 0.45 | 0.5 
Laa -0.3% | -0.3% | -0.2% | -0.2% | -0.1% | 0.0% | 0.1% | -0.2% | 0.0% | 0.1% | 0.3% 
Lab -1.2% | -1.2% | -1.2% | -1.1% | -1.0% | -0.9% | -0.8% | -1.1% | -0.9% | -0.7% | -0.6% 
Lac -1.2% | -1.2% | -1.2% | -1.1% | -1.0% | -0.9% | -0.8% | -1.1% | -0.9% | -0.7% | -0.6% 
Laf -3.2% | -3.2% | -3.2% | -3.1% | -3.1% | -3.0% | -2.9% | -3.1% | -3.0% | -2.9% | -2.7% 
Lf -3.1% | -3.1% | -3.1% | -3.0% | -2.9% | -2.8% | -2.7% | -2.9% | -2.8% | -2.6% | -2.4% 
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2.3.2 Simulated Results 
The SPSM behavior is simulated during a steady state condition. In this section, 
the simulation results of open circuit voltage waveforms and loaded current waveforms 


with DC field current source are presented and the spectrums with different eccentricity 


levels are compared. 


2.3.2.1 Open Circuit 


Stator current in time domain f=0 EccS=0 EccD=0 


Amplitude/V 


Figure 11 Voltage waveforms for the open circuit case 


The three phase voltages V,,V,,and V. at the steady state for the open circuit 
configuration are shown in Figure 11. Since the circuit is open, the terminal voltages are 
the same as the back-EMF of the circuit. Note that the voltage waveforms are more like a 


staircase than a sinusoidal shape, because the shape of back-EMF is determined by the 
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winding configuration. The staircase waveform represents a full pitch setup of the 
winding. 

Spectrums of the voltages with different dynamic eccentricities are shown in 
Figure 12. The simulation results from the WFA-based model with static eccentricity 
levels 0, and 50% are shown in Figure 12 (a), and (b), respectively, and the results for the 
FEM model are shown in Figure 12 (c) and (d)respectively. Similarly, the results for 
dynamic eccentricity are shown in Figure 13. For most of the harmonics listed in Table 
3, implementing SE and DE increases the induced voltage harmonics. The 19" harmonics 
obtained with FEM decreases while eccentricity is applied, which may due to the low 
resolution caused by short simulation period in FEM. The difference between the 
signatures of SE and DE is that 3“ harmonics is generated in the case of SE. It is because 


that the air-gap variation for each winding is not identical during one rotate of the rotor. 


Table 3 Relative increase of voltage harmonics due to 50% SE and DE — no load 
Order of harmonics 3rd 5th 7th 11th 13th 17th 19th 


e tee sas) dasha 
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Figure 12 Open circuit phase voltage spectrums - SE 
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Figure 13 Open circuit phase voltage spectrums - DE 
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2.3.2.2 Loaded with Resistors 

To simulate the loaded case, the synchronous generator model is connected with 
50ohm resistors in a star configuration. The field winding is fed with 4.8A DC current. 
Figure 14 and Figure 15 show that the current waveforms simulated from both models are 
not perfect sinusoidal signals. The distortion caused by the interaction of the magnetic 
fields produces harmonic fluxes that move relative to the stator, and induce corresponding 
current harmonics in the stationary stator windings. The simulated torque waveform of the 
WFA model covers the envelope of the FEM result. The torque waveform of the FEM 
model has higher frequency components than the WFA model, because the slot effect is 
ignored in the WFA model. 

The percentage changes of the harmonics for SE and DE with loaded case are 
shown in Table 4. Since the modeled SPSM is a three phase system which neutral is not 
connected to ground, the 3'* harmonic and its multiples are canceled out in the stator 
windings. The 5*' and 7" harmonics decrease while the eccentricity level increases in both 
models, whereas the 11", 13", 17", and 19" harmonics increase in both models. Either/or 
both signature analyses can be utilized to detect eccentricity problems, since the stator and 


rotor current signatures of synchronous machines have changed. 
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Compare FEM results and MWFA results at steady state 


Figure 14 Currents signal in time domain 
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Figure 15 Electromagnetic torque in time domain 
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Table 4 Percentage change of current spectrum due to SE and DE -loaded 


Order of harmonics 3rd 5th 7th 11th 13th 17th 19th 
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CHAPTER III 

FLEXIBLE ROTOR WITH ROLLING ELEMENT BEARING FAULTS 

In Chapter III, the methods of modeling faults of ball bearing are described and a 
nonlinear bearing, simplified fault model based on Hertzian contact theory is developed. 
A flexible rotor with built-on components is modeled using FEM. The supporting system 
between the housing and the base is modeled as linear springs and dampers. The excitation 
force caused by the bearing defects is treated as nonlinear. Numerical simulations are 
conducted to investigate the effects of two types of bearing faults, including outer raceway 


and inner raceway faults. 


3.1 Rolling Element Bearing Fault Model 

3.1.1 Root Causes of Rolling Element Bearing Failures 

The causes of bearing failures are perhaps as numerous as the bearing installations 
themselves. Tallian [67] and Johnson [68] listed the many modes of failure for rolling 
element bearings, all of which are preventable except contact fatigue. The cause or mode 
of failure appears as changes or irregularities in the rolling surfaces of the bearing. 

Contact fatigue is a type of a surface defect or damage that is inevitably related to 
the operational wear of rolling element bearings. Spalling, pitting, or flaking of the 
metallic particles from a bearing’s rolling surfaces, namely the outer raceway, inner 
raceway, and rolling elements, characterize contact fatigue. 

Faulty installation, including excessive preloading in either radial or axial 


directions, also causes bearing failure. Misalignment, loose fits, or damage occur because 
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of the excessive force used in mounting the bearing components. A uniformly wide wear 
track at the rotating raceway that extends over the entire circumference characterizes the 
misalignments. 

Excessive wear results from dirt and foreign particles entering the bearing because 
of inadequate lubrication. Severe wear changes the raceway profile and alters the rolling 
element profile, increasing the bearing clearance. An increase in rolling friction leads to 
high levels of slip and skidding, and eventual breakdown. 

Overheating, the ring and ball materials can be annealed by a temperature in excess 
of 400°C can anneal. The resulting loss in hardness reduces the bearing capacity and leads 
to premature failure. Discoloration from gold to blue of the rings, balls, and cages 
characterizes overheating. 

3.1.2 Bearing Model 

Studies of localized faults have used two main types of fault models. One method 
simulates the defect as a signal function of impulsive train for the modelled system. This 
impulsive train model treats the bearing as two sets of linear springs and dampers on the 
horizontal and vertical directions. This type of model requires less information and less 
computation effort, because of the linearly modeled bearing. The other method models the 
bearing defect as a surface bump or a dent with length, width, and depth. A high degrees 
of freedom (DOF) model is usually used for the bearing because the passing events of 
each rolling element,. The DOF is related to the total number of rolling elements. This 
type of model provides more possibilities for parametric studies, although the computation 


cost is greater. 
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The next sections describe a multiple DOF, simplified fault model of a ball bearing 
and its application to investigate the vibration of a faulty bearing running at a given speed. 
The assumptions are as follows: 

e All rolling elements, and raceways have motions in the plane of the bearing only. 

e The elastic deformations of the contacted surfaces are small. 

e The mass of inner ring and outer ring is ignored. 

e Deformations in the contact occur according to the Hertzian theory of elasticity. 

e The cage is rigid. Therefore, all rolling elements have the same angular speed. 

e The axial load is large enough, i.e., there is no skidding on the contact surface 
between the rolling elements and the raceways. 

e The flexible housing is modeled by the springs and dampers that connect nodes of 
shell elements to the rigid base. 

A 29 DOF nonlinear and time varying rotor bearing model is built in [67], where 
the rotor and housing are assumed to be rigid. In this research, the rotor and housing is 
treated as flexible and modeled using FEM. Lagrange method is used to derive the 


equations of motion for the ball bearing: 


cas) OT = 520 - (f) 3.1) 
dt\d{p}/ dtp} dtp} l 
where T is kinetic energy, V is potential energy, p is the vector of generalized DOF 
coordinate, and f is the vector of generalized force. 


To calculate the potential energy of the contact deformation, the displacements of 


ball, inner race, and outer race as shown in Figure 16 have to be computed first. The 
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contact force occurs only if there is a compression in the contact surface. Thus, the inner 


and outer contact deformation for the jth ball can be expressed as: 


asp en ifr+pp—-xX >0 
à 0, otherwise 
(3.2) 


Oy = (o t py) — Ro, if (pj +P») - Ro > 0 
0, otherwise 
where r is the radius of the inner race, p; is the radial position of the jth ball, p; is the 


radius of the ball, and the subscriptions i and o are the inner race and the outer race, 


respectively. 


Figure 16 The contact deformation between jth ball and the races 
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The offsets of the inner race and the outer race are given by: 
7, = xii t yij 
To = Xoi  yij (3. 3) 


p; = pjcosOji + pjsin8jj 


Since y; = r; — ri - p,, by substituting the expressions of r5,7,p, into the 
definition, the distance between the inner race center and the ball position can be obtained 


as: 


1 
x; = ((xo — xi + pjcos8j) + (Yo — yi + pjsin&j) Y (3. 4) 


Following the method in [59], the total potential energy of the bearing is given by: 


Np 

1 1 

— 2.5 2.5 T 

¥ = 75D (el 825) € zu Ku (3.5) 
J= 


where N, is the number of rolling elements in the bearing, u is the nodal displacement 
vector [Xi Yi Xo Yoļ",K is the stiffness matrix of the inner ring, the outer ring and 
the supporting base, ô; j; and 6, ; are the inner and outer deformation of the jth ball, k; and 


k, are the stiffness at the contacts of inner and outer races. K is calculated by Hertzian 


Contact Theory as: 
0 0 0 0 
_{0 0 0 0 
K-l|g 9 kns 0 (3. 6) 
00 0 ky 


The kinetic energy of jth rolling element is given by: 
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1 are > hee 1,, 


The first items in equation (3.7) are from the rotation energy of individual rolling 
elements with respect to the inertial frame. The second term represents the kinetic energy 
of rotating around the rotor. For simplicity, the ball spinning energy term is ignored. 


Substituting (3. 3) into (3. 7), Tp ; can be expressed as: 


1 H ^ Bj x 5 v . 
Taj = 5M (6; + p20? t x2 + yd + 2X, (pjcos0; — p;O;sin0;) 


(3. 8) 
+ 297 (pj6jcos0; + p;sin0;)) 
Because there is no mass at the nodes of the housing element and the shaft element, 


only the kinetic energy of the rolling elements must be considered. Thus, the total kinetic 


energy is obtained as: 


Np 
P X Pay (3.9) 
j=l 


The generalized coordinates are pj, where j = 1,2, .., Ny, xi, yj, xo and yo. 


Applying Lagrange method, the equations of motion are given by: 


mp; — pjc? + my (X,cos6; + y,sin8;) + koô? 


Ls Dj + (xo — xi)cos0; + (y, — yi)sin 0j (3. 10) 
E kiójj Ie NEC DE: deae MM tit (d CNET 0 
; X 
Np 2 
p;cos0; + Xo — xj 
ki i$ ————— — = Qsx (3. 11) 
j=1 Xj 
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pjsin8j + y, 
ki gps OTE - pcs 

jzi M 

Np 
NpMpXo + > mpp;cos0; + kyyx, 

j=1 

EX d TX» — Xi - 0, 
x 

Np 
Nympy, + »3 mppj;sin0; + kyyyo 

j=1 


I 2X, ax spem RA. 


the contact forces. The forces are the same but in the opposite direction. 


3.1.3 Modeling the Bearing Fault 


not all of the rolling elements roll over a bearing defect at any given time. 
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(512) 


(3. 13) 


(3. 14) 


where Q, and Q; are the virtual load applied on the node of the shaft element and the 


housing element, respectively. By comparing (3. 11)-(3. 14), the terms with 6; ; are clearly 


The bearing faults emulated in the experiments are small oval pits or holes 
manufactured into the raceway surfaces at specific locations. Adding time varying, 
position dependent raceway diameters is required to incorporate the localized effects of 
these bearing faults into the dynamic model. In the case of both inner and outer raceway 


faults, the bearing defect encompasses a small region of the raceway surface only; thus, 


The shapes of the localized defects have been discussed in many papers. The 
objective of this research is to develop a rotor bearing model which can generate the 
vibration signatures of the bearing faults so they can be used for the SPSM model. The 
relationship between the vibration signatures and the electrical signatures could be 
investigated. Therefore, modeling the shape of the defects is beyond the scope of this 


research. A simplified fault model is used in this dissertation instead. 


Figure 17 Diagram of outer race fault 
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A defect on the surface of the outer race at an angle $ from the horizontal axis, is 
shown in Figure 17. The defect formed by an oval pit takes over an angle length of Af 
on the surface from hq to dg + Ada. When the jth rolling element passes through this 
area, the contact deformation is recovered by the extent of the pit depth. Therefore, the 


deformation function (3. 2), can be replaced by: 


0, otherwise 


Once the rolling element passes the defect area, the deformation function is reset to (3. 2). 
In the case of an inner race defect, the defect rotates with the angular speed of the 
shaft w,. If the defect angle [$, $4] coincides with one of the balls, the deformation of 


that ball is written as: 


rtpp-Xj-D,  ifrtpp-x; -Di 20 


gines 0, otherwise (0:49) 
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3.2 Modeling the Rotor Bearing System 
In this section, the method of modeling and assembling a rotating rotor with 
flexible housing and nonlinear bearing forces is described. The rotor with built-on 
components is modeled linearly. The supporting bearing system is modeled nonlinearly 
as discussed in section 1. 


3.2.1 Flexible Rotor Modeling 


Bearing #1 


Housing and stato[ core 


Spring and damper Spring and damper 


FEREWIEPE TB P 
Figure 18 Diagram of the rotor-bearing system 


To investigate the vibrational motion, FEM is used to model the rotor bearing 


system. The shaft is assumed to be axial symmetric and modeled using a 4 DOF beam 
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element based on Timoshenko beam theory, considering shear deformation, rotary inertia, 
and gyroscopic effect. A 2 DOF shell element is used to model the non-rotating housing 
connected to the rigid base by springs and dampers. The rotor and stator of the motor are 
modeled as lumped masses located on the shaft and housing nodes. The governing 
equations of this rotor bearing system are expressed as: 


Mü + (C + QG)ù + Ku = F,(t) + FK, (u, ù, t) (3. 17) 


T 
Xi Yi Oxi Oyi 


where the displacement vector u = 
Xhj Ynj `“ Pk 


, Xi, yi Oxi , and 


Oy; represent the four DOFs of shaft node i, x; j, yn; represents the 2 DOF of housing 
node j, and p, is the radial offset of the rolling element k in the bearing, M is the mass 
matrix, G is the gyroscopic matrix, C is the damping matrix and K is the stiffness matrix; 
M,G,C, and K are linear terms. The nonlinear terms are moved to the right of (3. 17) and 
included in the nonlinear force F,(u, ù, t). E,(t) is the imbalance force caused by the 
imbalance mass, and Q is the angular velocity of the rotor. 
3.2.2 Assembling Nonlinear Bearing Fault Model 

The entire rotor-bearing system is separated into linear and nonlinear subsystems. 
First, the three subsystems (the flexible shaft with rotor disk, the nonlinear bearing 
supporting system, and the flexible housing with an attached stator part of the modelled 
motor) are modeled separately and then they are assembled. The four steps required to 
assemble the rotor bearing system are as follows: 

1) Generate the linear shaft and the built-on rotor disk matrices of each element using 

1D Timoshenko beam theory. 
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2) Assemble the linear rotor including the shaft and all built-on components and add 
it to the global matrices. The size of each matrix of the element is 4x4. 

3) Assemble the linear housing part and add it to the global matrices. The size of 
element matrix is 2x2. 

4) Add the nonlinear bearing system to the global matrices. 

Assuming the bearing is connected to the shaft at node i and to the housing at node 

j,the nonlinear bearing forces calculated from (3. 11)-(3. 14) are added to force vector 


at node i and j as: 


EC (3.18) 


where [K;]£h., is the shaft stiff coefficients of node i, [KH J is the housing stiff 
2x2 


coefficients of node j, and Fprg is the nonlinear bearing force. 


3.3 Simulation Results 


To validate the model, three cases are simulated: healthy machine with unbalanced 


mass only (Healthy Machine), unhealthy machine with bearings of the outer race defects 


(Outer Race Fault), and unhealthy machine with bearings of the inner race defects (Inner 


Race Fault). Experimental data from the Bearing Center at Case Western Reserve 


University are used to construct the test dataset. All parameters used in this research are 


obtained from [67]. The test rig is built by a Reliance Electric IQ Prealert 2-horse power, 


three phase induction motor. The original motor bearings are replaced by defected 
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bearings; the defects are added to the inner and outer raceways by using the electro- 
discharge machining (EDM) process. The defect size and depth are 0.007 inches and 0.011 
inches for the outer race defects and the inner race defects, respectively. Motor housing 
vibration data are measured at the 12 o'clock position of both the drive end and the fan 
end. The geometry of the bearings and the characteristic defect frequencies provided by 


the OEM are reported in Table 5 and Table 6, respectively. 


Table 5 Bearing geometry (in inches) [67] 


Bearing Inner Outside Thickness Ball Pitch 
Diameter Diameter Diameter Diameter 

Bearing 1 0.9843 2.0472 0.5906 0.3126 1.537 

Bearing 2 0.6693 1.5748 0.4724 0.2656 1.122 


Table 6 Bearing defect frequencies, where fm is the shaft rotational frequency [67] 


Bearing Inner Race | Outer Race | Cage Rolling 
Element 


Bearing 1  |5.4152f,  |3.5848f,  |03983f, | 4.7135fm 


Bearing2 | 4.9469f,, | 3.0530f,  |03817f,  |3.9874f,, 


3.3.1 Validation of the Rotor Bearing Model 

The orbits/airgap variations are plotted in Figure 19, Figure 20, and Figure 21. In 
a health machine with unbalanced mass, the centrifugal force caused by the mass would 
pull the shaft off its rotating center, which resulted in an elliptical shaft orbit. In the case 
of faulty bearings, the vibration caused by the bearing defects is not synchronized with the 


rotating speed. Thus, shaft orbits move chaotically, not in a perfect ellipse anymore. 
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Figure 19 Healthy Machine: Orbit Plot Simulation 
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Figure 20 Outer Race Fault: Orbit Plot Simulation 
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Figure 21 Inner Race Fault: Orbit Plot Simulation 


The housing vibrations are compared in both the time domain and the frequency 
domain. Waveforms of the housing vibration in the experimental data and the simulated 
data are shown in Figure 22-Figure 25. In the time domain, both waveforms illustrate the 


periodically beating phenomenon caused by the defect correctly. 
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Housing vibration with outer race defect - experiment 
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Figure 22 Waveform of experimental data with outer race defect 


Housing vibration with outer race defect - simulation 
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Figure 23 Waveform of simulated result with outer race fault 


61 


Housing vibration with inner race defect - experiment 
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Figure 24 Waveform of experiment data with inner race defect 
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Figure 25 Waveform of simulated result with inner race fault 
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Figure 26-Figure 29 show the spectral of the housing vibration. Because of the 
amplitude modulation, the slot between harmonics in the spectrum represents the defect 
characteristic frequency. The fault frequencies are about 110Hz for the outer race defect 
and 160Hz for the inner race defect, which match the findings on the spectrums. There are 
also qualitative similarities between the spectral of the real and simulated data. The 
bifurcation is about 3kHz in the real data spectrum and about 2.8kHz in the simulated data 
spectrum. The difference is acceptable considering the parameters’ error between the 


model and the real system. 
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Figure 26 Experiment spectrum of housing vibration under outer race fault 
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Housing vibration spectrum with outer race defect - simulation 
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Figure 27 Simulation spectrum of housing vibration under outer race fault 
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Figure 28 Experiment spectrum of housing vibration under inner race fault 
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Housing vibration spectrum with inner race defect - simulation 
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Figure 29 Simulation spectrum of housing vibration under inner race fault 


The simulation results are compared to the experiment results in Figure 30 and 
Figure 32. Both the peak values and the excitation frequencies are well-matched between 
the simulation and the experiment measurements. Vibration RMS of the waveforms are 
calculated and the results are summarized in Table 7. The acceptable errors are 2.2% for 


the outer race fault and 6.796 for the inner race fault. 
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Figure 30 Housing vibration with outer race fault 
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Figure 31 Housing vibration with inner race fault(zoomed) 
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Figure 32 Housing vibration with inner race fault 


Table 7 Vibration RMS 


Simulation | Measurement Error 


2.522 2.364 
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0.035 


3.3.2 Simulation Results of the Coupled Models 

With the similarities of the simulated data and the experiment data, the developed 
model is considered as valid and credible. The rotor bearing system model is coupled with 
the synchronous machine model to investigate the electrical signature changes. The 
flowchart of the coupled simulation is shown in Figure 33. At time kAt, the mechanical 
rotation speed wm, the rotation angle Om, and rotor center offset (x, y) are calculated 
based on the accelerations and velocities calculated at the last time step. The inductances 
are retrieved according to the position information. Electromagnetic torque T,, terminal 
voltage vector V, and phase current vector I are calculated by solving the electrical model 
equations. A trial solution for iteration step i is obtained by applying the boundary 
conditions to the mechanical model. If the solution error is less than the tolerance error, 
the variables and the time step are updated. Otherwise, the equilibrium iteration is re-run 
with the latest known solution. The process is stopped when the maximum time step tmax 


is reached. 
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Figure 33 Flowchart of the coupled simulation 


The current spectrums of the coupled system simulation with outer race fault and 
inner race fault are shown in Figure 34 and Figure 35, respectively. The existence of any 
bearing fault will result in a periodic airgap variation. The airgap changes will induce 
sideband components modulated on the harmonics at frequencies kf, + mfz, where k and 


m are integer numbers, f, is the vibration excitation caused by the mechanical fault, and 


fs is the line frequency. In the simulation, f. = 110Hz for the outer race fault and f; = 
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60Hz for the inner race fault. The peaks at the series of frequencies, which are clearly 
shown in the spectrums, reveal the interaction of the mechanical vibration and the 
electrical signature. According to the findings, the bearing faults’ information is included 
in the electrical signals, and therefore the electrical signals can be used as a source to 


perform the fault diagnosis for the mechanical faults. 
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Figure 34 Coupled simulation results with outer race faults 
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Current spectrum with inner race defect - simulation 
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Figure 35 Coupled simulation results with inner race faults 


To investigate the interaction between the severity of bearing defects and the 
electrical signatures, the amplitudes of frequency components along with the maximum 
deflection during both OR and IR are plotted in Figure 36 and Figure 37. Clearly, all fault 
related harmonics such as 169Hz, 369Hz in OR and 102Hz, 222Hz in IR increase when 
the vibration induced by the bearing fault goes up. However, unlike the case in SE and DE 
discussed in Chapter II, no observable pattern can be found on odd harmonics at 780Hz 


and 900Hz. 
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Figure 36 Frequency components changes by vibration level - OR 
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Figure 37 Frequency components changes by vibration level - IR 
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CHAPTER IV 
ESA AND EXPERIMENT RESULTS 
In this chapter, first the in-situ bearing deterioration experiment setup will be 
simply illustrated. Then we will investigate the nature of the generalized roughness 
bearing defect and its physical link to current signature. The power quality variation will 
be illustrated to show how the non-stationarity of power input could affect the ESA of 
subtle current signature of generalized roughness bearing defect. The discussion on the 
nonstationary, locally and piecewise stationary signal lead to applying principal 
component analysis (PCA) to the experimental data. The data processing results are shown 


in Section 3. 
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4.1 Experiment Setup 

4.1.1 In-situ Bearing Deterioration 

Modern AC converters have a high du/dt combined with a high switching 
frequency, i.e. the sum of the 3 phase voltages is no longer zero. This nonzero and high- 
frequency common-mode voltage can be considered the root cause for different kinds of 
bearing currents. Basically, the damage is always caused by partial discharge, which is 
called Electrical Discharge Machining (EDM). This film of lubricant serves as a dielectric 
and allows the bearing to behave as a capacitor. The capacitor created within a bearing 
between the inner /outer race and the ball is shown in Figure 38. This capacitor is charged 
by the bearing currents. When the voltage level is high enough, it is discharged by a short 


circuit. Eventually, the metal is eroded by the periodic discharging. 


Figure 38 Capacitor model between the bearing housing and the rolling element 


In [60], a shaft current injection method was proposed to emulate the generalized 
roughness bearing defect, which is a more practical anomaly in the field, and avoid the 


disassembly of the electrical machines and operational disruption. Therefore, the possibly 
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misleading vibration and electric signatures induced by reinstalling the bearing in the 
offline defect creation method are eliminated. 


4.1.2 Experiment Setup 
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Figure 39 The experiment setup of the in-situ bearing deterioration process 


The setup for the experiments performed for airgap variation induced by bearing 
generalized roughness is shown in Figure 39. Motor terminal quantities are sensed through 


current/voltage transducers and fed into a customized data acquisition box, and the 
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vibration signals are measured by an NI DAQ device. The two Marathon Electric 
synchronous machines 281PDL1722 with 8.8K W rated power, 240V/3 Phase are used in 
the test-bed as shown in Figure 40. One is used in motor mode as the primary drive and 


the other one is used in generator mode loaded by a water heater. 


(DANGER 


"Ay 


Figure 40 The testbed 


For current injection, an external voltage source is applied to the shaft of the motor 
via a carbon brush. Because of convenience and availability, a single-phase, 120V 60Hz 
AC power adjusted by a transformer is used to provide the controllable shaft voltage. 
Replacing the NDE bearing with a ceramic one allows a current to flow from the shaft 
through the drive end bearing and return to the power negative from the stator frame. The 


deterioration process is accelerated by removing part of the lubrication grease in the 
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bearing. An aluminum disk is mounted onto the shaft to provide a smooth contact surface 
for the brush. 

The vibration signals on the top of the bearing housing are collected to evaluate 
the deterioration level of the bearing. The accelerometer 352C33 with an amplitude range 
from -50g to 50g, a sensitivity of 100mV /G , and a frequency range from 0.2 — 5000Hz 
is used to measure the vibration, and the. DAQ system NI USB-9234 is used to collect 
the vibration data. The latter system, a USB-based four-channel C Series dynamic signal 
acquisition module that makes high-accuracy audio frequency measurements from 
sensors, is ideal for a wide variety of portable applications, such as industrial machine 
condition monitoring, in-vehicle noise, vibration, and harshness testing. The device and 


its chassis CDAQ-9174 is shown in Figure 41. 


Figure 41 Vibration DAQ system 
78 


Because the electrical signals cannot be measured directly by the data acquisition 
system, the phase currents and the line voltages have to be converted to low measurable 
voltage signals through current transducer (CT) and potential transducer (PT), respectively 
as shown in Figure 39. The AD73360, a six-input channel analog front-end processor, is 
used in the electrical signals collection. It is particularly suitable for industrial power 
metering as each channel samples synchronously, i.e. in this setup, the voltage signals and 
the current signals of three channels are sampled simultaneously. The sampling rate is set 


to 8 kHz. 


4.2 Electrical Signature Analysis 

In Chapters II and IV, the physical link between the airgap variation and current 
signatures in a synchronous machine is investigated. Based on the simulation results, the 
airgap variation caused by the bearing defects is modulated on the fundamentals and the 
odd harmonics of the electrical signals. Thus, the bearing defects are detectable by using 
electrical signature analysis. 

In a real-world environment, an incipient bearing defect usually occurs in the form 
of a generalized roughness surface, but it is not possible to identify a specific frequency 
in the defect's signature. Instead, unpredictable broadband changes in vibration and stator 
current spectrum are induced, which can be detected by simply using Root Mean Square 
(RMS) analysis. The RMS vibration value provides a good indicator for the health level 
of the machine, as stated in ISO-10816. However, in the electrical signature, the 


components after modulation are small compared to the fundamental harmonics. 
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Therefore, appropriate preprocessing methods and high-resolution data processing 


algorithms are required to analyze the electrical signature. 


4.2.1 Signal Segmentation 

In real world, current and voltage signals in an electrical machine are 
nonstationary. When electrical machines are running, factors such as power supply quality 
and load disturbances cause deviations in the values of current and voltage signals. For 
example, the load-induced harmonics in the stator current spectrum would obscure certain 
fault signatures, if the load of an electrical machine fluctuates over time. The sensitivity 
of any proposed ESA based fault detection methods would be degraded in this situation. 
When the measured signal is nonstationary, more complicated signal processing 
techniques are required in a fault diagnosis system than in a system dealing with stationary 
signals. 

Because of the nonstationary nature of the electrical signals, methods considering 
whole-spectrum energy content have been proposed to detect generalized roughness 
bearing faults, e.g., in [61] non-bearing faults components were removed from the signal 
and the residual was used to seek the frequency bands with a high probability of the 
presence of a bearing fault. The spectral kurtosis of the signals investigated in [62] [63], 
has proved useful for identifying the bandwidth where the fault-related signatures were 
more likely to appear on localized bearing faults and vibration monitoring. None of these 
studies, however, accounted for the influence of power quality change, load change, speed 


change, and other environmental variables. 
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Figure 43 Voltage THD variations with time 
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The voltage amplitude change during the 240 hours of normal operation is shown 
in Figure 42. The difference between the peak and the valley is about 2% of the average 
value. Total harmonics distortion (THD) is another measurement for the quality of the 
power input as shown in Figure 43.The voltage THD, which varies 25% during the 
operation, has significant impact on the magnitude of the modulated signal. The main 
sources of these variations may come from static frequency converters, cycloconverters, 
load changes in the power grid, etc. In the presence of these nonstationary variables, the 
electrical signature are difficult to pick up by directly applying spectral analysis approach 
[64]; thus, the data must be segmented in order to isolate the influence of the nonstationary 
variables. 

4.2.2 Feature Extraction 

Feature extraction is the essential step of the diagnostic of machinery, which can 
be treated as a problem of pattern classification/recognition. In the field of machinery fault 
diagnosis, signatures in time domain and frequency domain are the two main sources that 
features can be extracted from vibration analysis. The features both in time and frequency 
domain have been applied successfully [73]. In electrical signature analysis, because the 
fault signals are modulated on the line frequency and higher order harmonics, the features 
are usually extracted in frequency domain for stationary signal. In this dissertation, 
multiscale analysis is applied on the frequency domain and the features in the time domain 


are extracted by combining zoomed FFT and PCA. 
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4.2.2.1 Zoom Fast Fourier Transform 

As mentioned, the characteristic bearing fault frequency is modulated on the 
fundamental and space harmonics of the electrical signatures. Traditional FFT methods 
are implemented to perform the spectrum analysis of the stator current. The magnitude of 
these fault frequencies, however, is subtle compared to other frequency components. 
Therefore, to identify the subtle changes in the electrical signatures, a high-resolution 
frequency analysis method is required. 

The zoom FFT (ZFFT) is a useful algorithm to zoom in the narrow frequency band 
centered on a selected frequency [63]. The frequency resolution of the FFT is related to 
the sampling time and the number samples. In real applications on DSP-based systems, 
the size of FFT is limited by the memory available on the boards. Compared to traditional 
FFT method, ZFFT can reduce the sampling rate to decrease the number of samples and 
keep the same frequency resolution. The ZFFT algorithm used in this dissertation is shown 


in Figure 44. 
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Figure 44 ZFFT algorithm 


Initially, the sampled signal x[n] is obtained after sampling the signal x(t) at the 


rate f;. Then by multiplying e Ff ^, the spectrum of the signal is shifted by f. in the 
frequency domains as shown in (4. 1). In this way, the sampled sequence x,[n] is centered 
on fe. The complex signal x, [n] is then filtered by a low-pass filter with bandwidth Br. 
After filtering, the frequency range becomes [f., fe + By]. To extract the frequencies in 
this band, a 2B; sampling rate is enough according to the Nyquist-Shannon theorem. By 


is =, the new sequence x, [nM] has 
f 


downsampling the sequence at an integer rate M — 


N/M samples, and the reduced frequency signal dimension order saves both computation 


time and memory for the given bandwidth: 
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.Anmm N 


f fe] (4.1) 
The frequency spectrum includes the amplitude spectrum and power spectrum. 
Both are utilized in fault diagnosis. The odd harmonics and inter-harmonics are generated 
by different root causes. Odd harmonics in the electrical machine are usually caused by 
the setup of the windings, which is also called space harmonics, and synchronous 
vibrations such as oil whirl and rotor unbalance. The root causes of inter-harmonics can 
be bearing defects, shaft crack, gear faults, etc. Therefore, the frequency range is divided 
into 24 bins. Applying ZFFT to each bin obtains the spectrum of the expanded bandwidth. 
The sum of inter-harmonics in the selected range is calculated as fault features as well as 


odd harmonics. 


Table 8 The selected 24 features (f is the line frequency) 


Feature 7 
No. 


Of) fo (f2f) (2f,3f) 3f (3f,4f) (4f,5f) 5f ..  (15f,16f 
Servis MU (f,2f) (2f,3f) (3f,4f)  (4f,5f) ( ) 


4.2.2.2 Principal Component Analysis 

In vibration monitoring, the amplitude of vibration signals increases when the 
severity of the defect in the machine components degrades [65]. However, as the defect 
severity increases, the overall vibration curve is composed of an increasing trend and local 
variations of smaller magnitudes. Because of the stochastic nature of the fatigue failure 


progression process, it can be difficult to apply fault detection method to obtain usable 
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indicators from a contending feature set [78]. As shown in [76], features that provides 
incorrect information reduce the accuracy of fault diagnosis. Therefore, an appropriate 
feature selection method is needed to select the features that describe the defect condition 
best and eventually lead to reliable fault classification, diagnosis/prognosis. In 
multivariate data analysis, PCA is commonly used to analyze the covariance structure of 
data. This approach was developed to reduce the dimensionality of the input features for 
supervised and unsupervised classification purposes. 

Originated in the fundamental work of Pearson [66], PCA is generalized as a 
classic technique in statistical data analysis, feature extraction and data compression. The 
basic idea of PCA is that the elements are correlated mutually. Redundant information 
exists in multiple elements, which makes compression possible. By finding a rotated 
orthogonal coordinate such that the elements of x in the new space are uncorrelated, the 
redundancy induced by the correlation can be removed. 

Figure 45 shows that with PCA the sample is centered and rotated to line up with 
the direction of highest variance. It is noted that on the new coordinates, the variance on 


Z2 is small so that it can be ignored and the number of features is reduced from two to one. 
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Figure 45 Coordinate transform 


Considering a random vector x with n elements, a sample x(1),...,x(N) is 
collected by N measurements of x. Let X be the sample matrix [x(1) x(2) ...x(N)]. 
Based on PCA algorithm, the feature selection method is illustrated below: 

1) For the covariance matrix C = Cov(X), find a singular value decomposition C = 
UTSU such that S is a diagonal matrix and the singular values are put in 
descending order on the diagonal. 

2) Project the original features onto a new coordinate with the transform matrix U as 
y = Ux. 

3) Take into account the first k components that explain more than 90% of the 
variance according to X£.s;/X1.,s; , where s; is the ith largest singular value. 


4) Use the projected features to classify. 
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It should be noted that PCA only valid when the following assumptions and limits 
are met [81]: 
° The problem is linear. For nonlinear problem, other methodologies must be 


used to project the data from nonlinear space to the linear space such as kernel 


method. 
. The probability distribution of the variables must be normal. 
° The principal components are orthogonal, which offers an intuitive 


simplification that makes PCA soluble with linear algebra decomposition 
methods. 
. Large variances have significant dynamics, which implies the data have high 
SNR. 
4.2.3 Proposed Method 
As mentioned, the airgap variations caused by mechanical faults are included in the 
electrical signatures. ESA therefore does not require additional sensors and can provide 
continuous monitoring of the machine. It is also a non-invasive method. The vibration 
caused by bearing defects is modulated on the fundamental and higher order harmonics. 
Thus, the fault information is hidden in the sidebands of odd harmonics. In this 
dissertation, a new bearing faults diagnosis method that combines the high resolution 
spectral analysis method ZFFT with the PCA multivariate monitoring approach is 


proposed. The flowchart is shown in Figure 46. 
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Figure 46 Flowchart of the proposed method 


In the proposed diagnosis method, voltage and current signals for three phases are 
sampled simultaneously for every given time window t. ZFFT is then employed on each 
sideband of the odd harmonics, and the amplitudes of each frequency band are summed. 
The resulting metrics matrix with the features is reported in Table 8. Loading of the 
machine is used to segment the metrics matrix, PCA is applied on the segmented data, and 


the first two principal components are selected for the diagnosis. 
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4.3 Experiment Results 
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Figure 47 Housing vibration RMS for bearing deterioration experiment 


The experiment consists of 3 phases: a normal operation phase with health 
bearings, and two stages with defected bearings of different deterioration levels. Each 
phase lasts about four days. A current is injected into the bearing cage between phases to 
damage the bearings. The deterioration period lasts about one day. According to ISO 
10816, vibration RMS is preferred to evaluate the health status of the rotary machines. 
Thus, the vibration data is used both as the reference and to validate the ESA results. 

The RMS value of housing vibration during the entire experiment is shown in 
Figure 47. Each data sample is collected in a 60 sec interval. The trending of the vibration 


indicator gives the transition of health condition of the bearing during its operational 
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lifetime. As the severity level of the bearing fault increases, the vibration also increases. 
Three phases are clearly isolated by the vibration indicator. The vibration level is raised 
by 100% compared to the health case after each deterioration phase. The ramping up trend 
shows the deterioration process resulting from the current injection. During the current 
injection phase, the vibration indicator is not as smooth as the phase when the injected 
current is removed. The reason is that when there is a shaft current, the vibration pattern 
is affected by the thermal effect and the unbalanced magnetic pull caused by the injected 
current, both of which cause a larger variation on the ramping up trend of the vibration 
indicator. 

The vibration spectrums of the three phases are shown in Figure 48. The sampling 
frequency of the vibration DAQ is 12800Hz, according to which the bandwidth of the 
spectrum is 6400Hz, indicates that the vibration over a broadband increases as the 
condition of the bearing worsens. Because the vibration information is modulated on the 
electrical signals, the higher the frequency, the lower the amplitude of the modulated 
signal. 0—1000Hz is the frequency range of interest for ESA. The zoomed view of the 
vibration spectrum in Figure 48 (b) shows the series of peaks at 30Hz, 60Hz, 90Hz, ..., 
etc., induced by the unbalanced rotor. The peaks show little change when the bearing 
condition deteriorates. The broadband changes of the spectrum are from 600Hz to 800Hz 


over the frequency band. 
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Figure 48 Vibration spectrums 


During the long operation time, the electrical signal is sensitive to the environment 
variables and easily polluted by the loop current from the grid. The accumulated spectrum 
energy of the three selected bands which vary in different patterns is shown in Figure 49. 
Clearly, the time series of the variables are stochastic processes. The distribution of 
spectrum energy for the fundamental harmonics is mostly closer to the normal distribution 
as shown in Figure 49(b). The distributions of the spectrum energy for the 1 1th harmonics 
and the band (360Hz, 420Hz) are skewed to the lower values, which indicates that the high 


peaks in the time series are more likely to be noises. 
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Figure 49 Time series and histograms of some features 
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The band (360Hz, 420Hz) is examined further to investigate the changes with 
bearing deterioration. Zooming in on the base line of the series indicates that the low 
values of the series are clearly separated into 3 stages as shown in Figure 51. According 
to the histogram plot of the variables in Figure 49, 80% of the points are in the range of 
400 to 600, which are useful data for ESA. Moreover, if the whole dataset is used, the 
distribution of the variables eliminates the assumption that the probability distribution of 
x; must be Gaussian. Thus, data segmentation should be applied to the dataset prior to 


extracting the fault information. 
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Figure 50 Inter-harmonics of 360—420Hz without segmentation 
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Figure 51 Inter-harmonics of 360-420Hz without segmentation (zoomed) 


The correlation coefficients are calculated and summarized in Table 9. For 
simplicity, only part ofthe variables are listed. The correlation between variables is higher 
if the absolute value of the coefficient is larger. Note that all inter-harmonics to some 
extent are correlated in the same direction. The 11" harmonics is correlated in the opposite 
direction with the 7 and the 9" harmonics and has little correlation with the 3'* and the 
5 harmonics. Based on the analysis, the vibration information appears on multiple bands 
in the selections and the basic assumption holds by applying the data segmentation 
method. After pre-processing the data, it is preferable to use PCA to extract the features 


from the correlated variables. 
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Table 9 Correlation coefficients of inter-harmonics and harmonics 


[180- — [240- [300-  [360- [420- — [480- [540-  [600- 

240Hz  300Hz] 5th 360Hz]  420Hz] 7th 480Hz]  540Hz] 9th 600Hz] 660Hz] 11th 
3rd -0.03 -0.02 -0.04 -0.02 -0.01 i -0.01 -0.03 005 -003 0.02 -0.06 
[180-240Hz . . -0. 0.08 032 004 
[240-300Hz | 54 0.50 0.15 
Sth -0. -013  -0.19 -0. . -0. -0. : -0.16 0.07 
[300-360Hz 028 — -007 
[360-320Hz] — -0. 049 058 . . | . 47 ! 054 -0.05 
7th . -0.12 -0.20 -0.03 20883] 
[420-480Hz 038 -0.12 
[480-540Hz 035 — 001 
9th 0.30 
[540-600Hz 052 0.08 
[600-660Hz 002 032 O50  -016 028 O54  -003 038 035 030 05 0.04 
11th -0.06 004 015 007  -007  .O05[77 X081 -o12  — 001  -059 008 0.04 


After using PCA on the segmented data, principal components #1 and #2 are 
shown in Figure 52 and Figure 53, respectively. Comparing Figure 52, Figure 54, and 
Figure 55 shows the harmonics distortion induced from the power quality variations is 
represented by principal component #1. The changes of harmonics distortion come from 
the same source and have a variation larger than the inter-harmonics, where the bearing 
faults information are most likely hidden. 

The visually recognizable trend of principal component #2 is shown in Figure 53. 
The moving average operation is applied to smooth the change of the indicator. The trend 
shows that as the bearing deteriorates principal component #2 increases. The results prove 
that the electrical indicator created from the proposed method extracts the bearing faults 


from the electrical signals and shows a good agreement with the vibration indicator. 
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Figure 52 Principal component #1 
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Figure 53 Principal component #2 
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Figure 54 The 9th harmonics 
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CHAPTER V 
SUMMARY AND FUTURE WORK 
In this section, a summary of the research work discussed in earlier sections is 
presented. This dissertation was mainly concerned with developing a SPSM model which 
is good to analyze the relationship between the bearing defects and the electrical signatures. 


A new diagnosis technique based on the findings is also proposed. 


5.1 Summary of Research 

At the start of the dissertation, common failures of electrical machines are briefly 
overviewed. The causes of bearing anomalies and the detection techniques are discussed. 
An introduction of the development for bearing fault modeling is also presented. In the 
end, the research objectives and a roadmap of the proposed method are outlined. 

In the following chapter, a winding function approach based method was 
introduced to model the SPSM. In this model, the actual winding distribution is used, 
which includes the space harmonics caused by the winding configuration. The airgap 
variation can be taken into account in the proposed model, which can be used to simulate 
the machine under different bearing anomalies. All machine inductances were calculated 
as a function of the rotor position using the winding function method, based on the 
geometry of the machine and its winding configuration. A simple method of modeling the 
dynamic eccentricity was also presented. The electromagnetic torque equation was 
derived, and the machine differential equations were given which are solved using Matlab. 


The simulated machine waveforms and their spectrums were presented for the cases 


99 


during different level of dynamic eccentricity at the steady state, for open circuit mode 
and loaded mode. 

Then, to investigate the mechanism that how the bearing faults would affect the 
electrical signatures, a high fidelity rotor-bearing system with detailed nonlinear ball 
bearing model and a flexible finite element shaft model was built. Dynamic differential 
equation of the defect vibration model was numerically solved and the localized defect of 
the rolling element bearing was simulated by the ode solver in Matlab. The simulation 
result was theoretically analyzed and experimentally verified. The validity of the proposed 
defect model is thus verified. The simulation coupled both rotor bearing system and the 
electromagnetic circuit model of a SPSM was performed. It was found that the lateral 
vibration induced by bearing faults are modulated on the harmonics of the electrical 
signals, which could provide theoretical basis to the research on the potential condition 
monitoring and fault diagnosis algorithms of rolling element bearing in the future. 

Spalls over an extended area on the bearing race surface caused by inject shaft 
current are more practical than localized single point bearing fault. In Chapter IV, the 
experiment data for the general roughness bearing faults was analyzed and a new method 
to process the non-stationary signals was proposed by applying ZFFT and PCA. The 
results demonstrated the capability of using electrical signals to detect bearing defects in 
a SPSM. 

The main work of this research is to build an inter-coupled simulation model in 


bridging the gap between electromagnetic circuit modeling of the electrical system and 
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the dynamical modeling of the mechanical system. The following contributions are 
achieved through this research. 
e A multi-physics model combining a detailed bearing model in a synchronous 
generator is developed, which allows introducing bearing faults 
e SPSM voltages and currents under bearing defects induced air-gap variation and 
torsional vibration are analyzed 
e A new fault detection method is proposed, which is capable of picking up the 


bearing faults in a noisy environment. 


5.2 Suggestions for Future Work 
An inter-coupled model developed in this research bridges the gap between 
electromagnetic circuit modeling of the electrical system and the dynamical modeling of 
the mechanical system. The model has been validated with FEM and experimental data to 
be a viable tool to analyze rolling element bearing faults. The proposed fault detection 
algorithm for bearing faults has proven to be effective in an experimental study. However, 
there are still some improvements to be done to make the model and the approach more 
reliable and extendable to more general cases. 
e Damper winding and voltage fed field circuit could be taken into consideration, 
which may lead to the detection of failures through field current. 
e The interaction of the rotor and the stator due to the unbalanced magnetic pull 
(UMP) could be explored and modeled. 
e Shaft orbit profile could be obtained by Eddy current sensors through the 


experiment. This may provide a more direct way to validate the developed model. 
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e An ideal fault detection method should be load independent, so that it is more 
applicable for practical analysis. The influence of load changes on the proposed 


fault detection method should be examined in the future. 
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